
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% The impact of divorce laws on the equilibrium in the marriage market.
% Reynoso
% February 2024
%
% This file inputs outputs from the model estimated under MCD
% and the model simulated under UD and constructs marriage market
% and family life cycle outcomes.
%  
% Data: PSID 1968-1992
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

%------------------------ Preliminaries ----------------------------------%

%--- Working directory
cd ([replication_location,'\Outputs']);

%--- Path for estimates
estimates_dir = [replication_location,'\Outputs\MCD\'];

%--- Path to store
output_dir = [replication_location,'\Outputs'];

%--- Load outputs
load([estimates_dir, 'outputs_mcd']);
params=X;
load([output_dir, '\MCD\moments_mc']); 
load([output_dir, '\Divorce_settlement\moments_ud_DivorceSettlement']); 
load([output_dir, '\MCD\U_mc']); 
load([output_dir, '\Divorce_settlement\U_ud_DivorceSettlement']); 
load([output_dir, '\Divorce_settlement\equil_ud_DivorceSettlement']); 

%------------------------ UD with Divorce Settlement outcomes ------------%

%--- Organize moments
model_moments_bymkt_mc = reshape(model_moments_mc(1:168), size(model_moments_bymkt_ud,1),size(model_moments_bymkt_ud,2)); 
mf_f_mc = model_moments_bymkt_mc(1:9,:); 
mf_f_ud = model_moments_bymkt_ud(1:9,:); 
sahw_mc = model_moments_bymkt_mc(25:33,:); 
divorce_mc = model_moments_bymkt_mc(34:42,:); 
sahw_ud = model_moments_bymkt_ud(25:33,:);
divorce_ud = model_moments_bymkt_ud(34:42,:); 
sorting_f_mc = zeros(n_g+1,4);
sorting_f_ud = zeros(n_g+1,4);
for mkt=1:4
sorting_f_mc(:,mkt)=[(mf_f_mc(1,mkt)+mf_f_mc(5,mkt)+mf_f_mc(9,mkt))/sum(mf_f_mc(:,mkt)) ;
    mf_f_mc(:,mkt)./sum(mf_f_mc(:,mkt))];
end
for mkt=1:4
sorting_f_ud(:,mkt)=[(mf_f_ud(1,mkt)+mf_f_ud(5,mkt)+mf_f_ud(9,mkt))/sum(mf_f_ud(:,mkt)) ; 
    mf_f_ud(:,mkt)./sum(mf_f_ud(:,mkt))];
end

%--- Assortative matching measures
fraction_married_mc = sum(mf_f_mc);
hsm=(mf_f_mc(1,:) + mf_f_mc(4,:) +mf_f_mc(7,:))./fraction_married_mc;
scm=(mf_f_mc(2,:) + mf_f_mc(5,:) +mf_f_mc(8,:))./fraction_married_mc;
cpm=(mf_f_mc(3,:) + mf_f_mc(6,:) +mf_f_mc(9,:))./fraction_married_mc;
married_men_byeduc_bymkt_mc = [hsm;scm;cpm];
hsf=(mf_f_mc(1,:) + mf_f_mc(2,:) +mf_f_mc(3,:))./fraction_married_mc;
scf=(mf_f_mc(4,:) + mf_f_mc(5,:) +mf_f_mc(6,:))./fraction_married_mc;
cpf=(mf_f_mc(7,:) + mf_f_mc(8,:) +mf_f_mc(9,:))./fraction_married_mc;
married_women_byeduc_bymkt_mc = [hsf;scf;cpf];

clear hsm scm cpm hsf scf cpf

fraction_married_ud = sum(mf_f_ud);
hsm=(mf_f_ud(1,:) + mf_f_ud(4,:) +mf_f_ud(7,:))./fraction_married_ud;
scm=(mf_f_ud(2,:) + mf_f_ud(5,:) +mf_f_ud(8,:))./fraction_married_ud;
cpm=(mf_f_ud(3,:) + mf_f_ud(6,:) +mf_f_ud(9,:))./fraction_married_ud;
married_men_byeduc_bymkt_ud = [hsm;scm;cpm];
hsf=(mf_f_ud(1,:) + mf_f_ud(2,:) +mf_f_ud(3,:))./fraction_married_ud;
scf=(mf_f_ud(4,:) + mf_f_ud(5,:) +mf_f_ud(6,:))./fraction_married_ud;
cpf=(mf_f_ud(7,:) + mf_f_ud(8,:) +mf_f_ud(9,:))./fraction_married_ud;
married_women_byeduc_bymkt_ud = [hsf;scf;cpf];

%//random matchings

%-MC
rm_mc_mkt1 = married_women_byeduc_bymkt_mc(:,1)*married_men_byeduc_bymkt_mc(:,1)';
rm_mc_mkt2 = married_women_byeduc_bymkt_mc(:,2)*married_men_byeduc_bymkt_mc(:,2)';
rm_mc_mkt3 = married_women_byeduc_bymkt_mc(:,3)*married_men_byeduc_bymkt_mc(:,3)';
rm_mc_mkt4 = married_women_byeduc_bymkt_mc(:,4)*married_men_byeduc_bymkt_mc(:,4)';

%-UD
rm_ud_mkt1 = married_women_byeduc_bymkt_ud(:,1)*married_men_byeduc_bymkt_ud(:,1)';
rm_ud_mkt2 = married_women_byeduc_bymkt_ud(:,2)*married_men_byeduc_bymkt_ud(:,2)';
rm_ud_mkt3 = married_women_byeduc_bymkt_ud(:,3)*married_men_byeduc_bymkt_ud(:,3)';
rm_ud_mkt4 = married_women_byeduc_bymkt_ud(:,4)*married_men_byeduc_bymkt_ud(:,4)';

%//actual matchings

%-MC
am_mc_mkt1 = reshape(sorting_f_mc(2:10,1),3,3);
am_mc_mkt2 = reshape(sorting_f_mc(2:10,2),3,3);
am_mc_mkt3 = reshape(sorting_f_mc(2:10,3),3,3);
am_mc_mkt4 = reshape(sorting_f_mc(2:10,4),3,3);

%-UD
am_ud_mkt1 = reshape(sorting_f_ud(2:10,1),3,3);
am_ud_mkt2 = reshape(sorting_f_ud(2:10,2),3,3);
am_ud_mkt3 = reshape(sorting_f_ud(2:10,3),3,3);
am_ud_mkt4 = reshape(sorting_f_ud(2:10,4),3,3);

%//difference
ggks_mc_mkt1 = am_mc_mkt1 - rm_mc_mkt1;
ggks_mc_mkt2 = am_mc_mkt2 - rm_mc_mkt2;
ggks_mc_mkt3 = am_mc_mkt3 - rm_mc_mkt3;
ggks_mc_mkt4 = am_mc_mkt4 - rm_mc_mkt4;
ggks_ud_mkt1 = am_ud_mkt1 - rm_ud_mkt1;
ggks_ud_mkt2 = am_ud_mkt2 - rm_ud_mkt2;
ggks_ud_mkt3 = am_ud_mkt3 - rm_ud_mkt3;
ggks_ud_mkt4 = am_ud_mkt4 - rm_ud_mkt4;
ggks_mc_bymkt = zeros(3,3,4);
ggks_mc_bymkt(:,:,1)=ggks_mc_mkt1;
ggks_mc_bymkt(:,:,2)=ggks_mc_mkt2;
ggks_mc_bymkt(:,:,3)=ggks_mc_mkt3;
ggks_mc_bymkt(:,:,4)=ggks_mc_mkt4;
ggks_ud_bymkt = zeros(3,3,4);
ggks_ud_bymkt(:,:,1)=ggks_ud_mkt1;
ggks_ud_bymkt(:,:,2)=ggks_ud_mkt2;
ggks_ud_bymkt(:,:,3)=ggks_ud_mkt3;
ggks_ud_bymkt(:,:,4)=ggks_ud_mkt4;

%//Ratio of traces
delta_bymkt_ud=[trace(am_ud_mkt1)/trace(rm_ud_mkt1); trace(am_ud_mkt2)/trace(rm_ud_mkt2); trace(am_ud_mkt3)/trace(rm_ud_mkt3); trace(am_ud_mkt4)/trace(rm_ud_mkt4)];

%--- Gains from marriage:

EVM_mc_bymkt = log(U_mc_tot_bymkt) ;
EVMf_mc_bymkt = [EVM_mc_bymkt([1 4 7],1,1) EVM_mc_bymkt([1 4 7],1,2) EVM_mc_bymkt([1 4 7],1,3) EVM_mc_bymkt([1 4 7],1,4)];
EVM_ud_bymkt = log(U_ud_tot_bymkt)  ;
EVMf_ud_bymkt = [EVM_ud_bymkt([1 4 7],1,1) EVM_ud_bymkt([1 4 7],1,2) EVM_ud_bymkt([1 4 7],1,3) EVM_ud_bymkt([1 4 7],1,4)];
GM_mc = [log(U_mc_tot_bymkt(:,1,:)) - U_mc_bymkt(:,2,:)  log(U_mc_tot_bymkt(:,2,:)) - U_mc_bymkt(:,4,:)];
GMf_mc = [GM_mc([1 4 7],1,1) GM_mc([1 4 7],1,2) GM_mc([1 4 7],1,3) GM_mc([1 4 7],1,4)];
GM_ud = [log(U_ud_tot_bymkt(:,1,:)) - U_ud_bymkt(:,2,:)  log(U_ud_tot_bymkt(:,2,:)) - U_ud_bymkt(:,4,:)];
GMf_ud = [GM_ud([1 4 7],1,1) GM_ud([1 4 7],1,2) GM_ud([1 4 7],1,3) GM_ud([1 4 7],1,4)];
%// Consumption equivalence
Pi_c_f = (1- exp((EVMf_ud_bymkt-EVMf_mc_bymkt)/9.15))*100;

%--- Save
mf_f_ud_ds =mf_f_ud;
sahw_ud_ds =sahw_ud;
divorce_ud_ds =divorce_ud;
GMf_ud_ds = GMf_ud;
Pi_c_f_ds = Pi_c_f;
delta_bymkt_ud_ds = delta_bymkt_ud;
ggks_mc_bymkt_ds = ggks_mc_bymkt;
ggks_ud_bymkt_ds = ggks_ud_bymkt;

save('Divorce_settlement\UD_effects_DivorceSettlement',  'GMf_ud_ds',  'Pi_c_f_ds',  'delta_bymkt_ud_ds', 'ggks_mc_bymkt_ds', 'ggks_ud_bymkt_ds',  'sahw_ud_ds',  'divorce_ud_ds' ,  'mf_f_ud_ds'); 

